## A Research Note on the Prevalence
## of Housing Eviction among Children
## Born in American Cities

## CODE PART B:
## Stan code to fit models
## This file just has character strings to pass to stan,
## which stan turns into samplers.

## Ian Lundberg and Louis Donnelly

## Department of Sociology, Office of Population Research,
## and Center for Research on Child Wellbeing
## Princeton University

## Code by Ian Lundberg (ilundberg@princeton.edu)
## Last updated 8 May 2018

library(rstan)

stan_code = "
data {
// Define variables in data

int<lower=0> N_city; // Number of cities
int<lower=0> N_person; // Number of persons
int<lower=0> N; // Number of person-periods

// Cluster IDs for 
int<lower=1> city[N];
int<lower=1> person[N];

// Outcomes
int<lower=0,upper=1> Y[N];

// Number of predictors
int<lower=0> K;

// Predictors
matrix[N,K] X;

}

parameters {
// Define parameters to estimate

// Grand coefficients
real alpha;
vector[K] beta;

// Variance parameters
real<lower=0> sigma_city;
real<lower=0> sigma_person;

// Coefficients at each level
real u_city[N_city];
real u_person[N_person];
}

model {
// Variance parameters
sigma_city ~ cauchy(0,1);
sigma_person ~ cauchy(0,1);

// Betas at each level
alpha ~ cauchy(-4.5,1);
beta ~ cauchy(0,1);
u_city ~ normal(0, sigma_city);
u_person ~ normal(0, sigma_person);

// Outcome model
for (i in 1:N) {
Y[i] ~ bernoulli_logit(alpha + u_city[city[i]] + u_person[person[i]] + X[i]*beta);
}
}"

stan_code_noRE = "
data {
// Define variables in data

int<lower=0> N; // Number of person-periods

// Outcomes
int<lower=0,upper=1> Y[N];

// Number of predictors
int<lower=0> K;

// Predictors
matrix[N,K] X;

}

parameters {
// Define parameters to estimate

// Grand coefficients
real alpha;
vector[K] beta;

}

model {

// Betas at each level
alpha ~ cauchy(-4.5,1);
beta ~ cauchy(0,1);

// Outcome model
for (i in 1:N) {
Y[i] ~ bernoulli_logit(alpha + X[i]*beta);
}
}"

stan_code_noCovs = "
data {
// Define variables in data

int<lower=0> N_city; // Number of cities
int<lower=0> N_person; // Number of persons
int<lower=0> N; // Number of person-periods

// Cluster IDs for 
int<lower=1> city[N];
int<lower=1> person[N];

// Outcomes
int<lower=0,upper=1> Y[N];

}

parameters {
// Define parameters to estimate

// Grand coefficients
real alpha;

// Variance parameters
real<lower=0> sigma_city;
real<lower=0> sigma_person;

// Coefficients at each level
real u_city[N_city];
real u_person[N_person];
}

model {
// Variance parameters
sigma_city ~ cauchy(0,1);
sigma_person ~ cauchy(0,1);

// Betas at each level
alpha ~ cauchy(-4.5,1);
u_city ~ normal(0, sigma_city);
u_person ~ normal(0, sigma_person);

// Outcome model
for (i in 1:N) {
Y[i] ~ bernoulli_logit(alpha + u_city[city[i]] + u_person[person[i]]);
}
}"

stan_code_noREorCovs = "
data {
// Define variables in data

int<lower=0> N; // Number of person-periods

// Outcomes
int<lower=0,upper=1> Y[N];

}

parameters {
// Define parameters to estimate

// Grand coefficients
real alpha;
}

model {

// Betas at each level
alpha ~ cauchy(-4.5,1);

// Outcome model
for (i in 1:N) {
Y[i] ~ bernoulli_logit(alpha);
}
}"
